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ABSTRACT 

Context. In the last few years, the so-called "Nice model" has got a significant importance in the study of the formation and evolution 
of the solar system. According to this model, the initial orbital configuration of the giant planets was much more compact than the 
one we observe today. 

Aims. We study the formation of the giant planets in connection with some parameters that describe the protoplanetary disk. The aim 
of this study is to establish the conditions that favor their simultaneous formation in line with the initial configuration proposed by 
the Nice model. We focus in the conditions that lead to the simultaneous formation of two massive cores, corresponding to Jupiter 
and Saturn, able to achieve the cross-over mass (where the mass of the envelope of the giant planet equals the mass of the core, and 
gaseous runway starts) while Uranus and Neptune have to be able to grow to their current masses. 

Methods. We compute the in situ planetary formation, employing the numerical code introduced in our previous work, for different 
density profiles of the protoplanetary disk. Planetesimal migration is taken into account and planetesimals are considered to follow a 
size distribution between r 1 ™ (free parameter) and rj'"" = 100 km. The core's growth is computed according to the oligarchic growth 
regime. 

Results. The simultaneous formation of the giant planets was successfully completed for several initial conditions of the disk. We find 
that for protoplanetary disks characterized by a power law (X oc r~ p ), smooth surface density profiles (p < 1.5) favor the simultaneous 
formation. However, for steep slopes (p ~ 2, as previously proposed by other authors) the simultaneous formation of the solar system 
giant planets is unlikely. 

Conclusions. The simultaneous formation of the giant planets -in the context of the Nice model- is favored by smooth surface density 
profiles. The formation time-scale is in agreement with the estimates of disk lifetimes if a significant mass of the solids accreted by 
the planets is contained in planetesimals with radii < 1 km. 

Key words. Planets and satellites: formation - Planet-disk interactions - Methods: numerical 



1. Introduction 

The initial configuration of the Nice model (Tsiganis et al., 2005; 
Gomes et al., 2005; Morbidelli et al., 2005) represents the or- 
bital configuration of the outer solar system after the gas of the 
primordial nebula dissipated. The model assumes that the giant 
planets were initially in nearly circular and coplanar orbits. This 
is compatible with the work of Thommes et al. (2008), who find 
that planetary system analogues to our solar system come out if 
the gas giants do not undergo significant migration during their 
formation and remain in nearly circular orbits. The Nice model 
proposes a compact initial configuration for the location of the 
giant planets; more precisely, the giant planet system is assumed 
to be in the range of ~ 5.5 AU to ~ 14 AU. The gas giants, 
Jupiter and Saturn, are supposed to be close to their mutual 2:1 
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mean motion resonance (MMR); Jupiter located at ~ 5.5 AU 
and Saturn between 8 - 8.5 AU. This is an important condition, 
required to avoid the migration of both planets during the gas 
disk lifetime (Masset & Snellgrove, 2001; Morbidelli & Crida, 
2007). Regarding the ice giants, they are assumed to be located 
at ~ 1 1 and ~ 14 AU. 

Another important key of the Nice model is the existence of 
a planetesimal disk beyond the orbit of the giant planets. The 
inner edge of the disk is proposed to be located ~ 16 AU, and 
the outer edge is fixed at ~ 30 AU; the total mass of the plan- 
etesimal disk being ~ 35 - 40 M e . The planetesimals conform- 
ing this disk gravitationally interact with the giant planets and 
cause the inward migration of Jupiter and the outward migration 
of Saturn to their current positions. In this process, Jupiter and 
Saturn cross their mutual 2:1 mean motion resonance, leading 
the system to undergo a phase of dynamical instability. During 
this phase, Uranus and Neptune migrate outward chaotically. 
According to the numerical simulations of Tsiganis et al. (2005), 
there is a 50 % probability that the two icy planets switch places 
in the process. 

The success of the Nice model relies on the fact that it is able 
to explain quantitatively the final orbits, eccentricities, and incli- 
nations of the giant planets (Tsiganis et al., 2005); the chaotic 
capture of Jupiter's Trojan asteroids (Morbidelli et al., 2005); the 



2 



O. M. Guilera et al.: Simultaneous formation of Solar System giant planets 



origin of the Late Heavy Bombardment (Gomes et al., 2005); 
the formation of the Kuiper belt (Levison et al., 2008) and the 
secular architecture of the outer solar system (Morbidelli et al., 
2009). Remarkably, the Nice model does not say anything about 
how the initial conditions proposed by itself could be achieved. 

More recently, Morbidelli & Crida (2007) -following the 
work of Masset & Snellgrove (2001)- showed that under certain 
parameters of the gas disk, Jupiter and Saturn did not migrate 
once both planets got locked in 2:3 MMR. Morbidelli & Crida 
(2007) showed that this kind of no-migrating (or slowly migrat- 
ing) evolution is possible if the mass ratio of the two planets is 
similar to the mass ratio of Jupiter and Saturn. Then, Morbidelli 
et al. (2007) extended this study to the four giant planets of the 
solar system. They locked Jupiter and Saturn in 2:3 MMR, and 
found that the inner ice giant could be trapped in the 2:3 or 3:4 
MMR with Saturn. Then, the outer ice giant could be trapped in 
the 3:4, 4:5, or 5:6 MMR with the inner ice giant. They showed 
that the resonant structure is preserved until the gas is completely 
dispersed. Then, only two configurations were dynamically sta- 
ble for several hundred million years. Finally, they showed that 
the existence of a planetesimal disk beyond the outer ice giant 
lead to a dynamical instability, and the essential ingredients of 
the Nice model are preserved. Batygin & Brown (2010) extended 
this work and found more possible multi-resonant initial condi- 
tions, several of which evolved -in the presence of a planetesimal 
disk beyond the ice giants- preserving the essential ingredients 
of the Nice model. While this initial conditions are not exactly 
the same to that proposed by the Nice model, all of them propose 
a compact orbital configuration for the giant planets of the solar 
system in a similar way to the Nice model. However, all these 
studies assumed that the giant planets were already formed, and 
did not address the question of how they formed. 

Desch (2007) was the first to investigate this topic. Adopting 
the in situ formation and perfect accretion for the giant plan- 
ets, he recalculated the "minimum-mass solar nebula" consid- 
ering the initial more compact configuration of the giant plan- 
ets together with the remnant planetesimal disk, as proposed 
by the Nice model. He then derived a much steeper minimum 
mass solar nebula surface density profile (Z oc a -2168 , where a 
is the distance to the Sun) than the classical one (E oc a~ 3 / 2 , 
Weidenschilling, 1977; Hayashi, 1981). He also estimated the 
growth time for the cores of the four giant planets. Assuming the 
oligarchic growth regime for the solid embryos, relative veloci- 
ties in equilibrium, and a single size population of planetesimals 
with radii of 100 m, he obtained growth times of about 0.5 Myr 
for Jupiter, 1 .5 - 2.0 My for Saturn, 5.5 - 6.0 Myr for Neptune, 
and 9.5 - 10.5 Myr for Uranus. However, he did not consider the 
presence of the also growing gas envelope of the planets. 

In a more realistic way, adopting the initial configuration 
for the location of the giant planets as proposed by the Nice 
model (Jupiter at 5.5 AU, Saturn at 8.3 AU, Neptune at 1 1 AU, 
and Uranus at 14 AU), and employing the surface density pro- 
file derived by Desch (2007), Benvenuto et al. (2009) calculated 
the in situ, isolated formation of the giant planets of the solar 
system. They used the same code already introduced in pre- 
vious works (Benvenuto & Brunini, 2005; Fortier et al., 2007; 
2009) but incorporating a size distribution for the planetesimals 
radii, considering they take values between 30 m - 100 km (nine 
species geometrically evenly spaced), with a mass distribution 
n(m) oc mT a , with a = 2.5. The choice of r^ m = 30 m as the 
minimum radius of the size distribution of planetesimals relied 
on a simple estimation for the planetesimal radius for which the 
migration time-scale of smaller planetesimals should be shorter 
than the protoplanet's planetesimal accretion time-scale, so the 



accretion of planetesimals with radii < 30 m should be not so ef- 
ficient. Imposing a value of 1 1 g crrT 2 for the solids surface den- 
sity at the location of Jupiter, they calculated the formation time 
for Jupiter (0.44 Myr), Saturn (1.4 Myr), Neptune (2.5 Myr), 
and Uranus (4.75 Myr). Furthermore, they found that the final 
core masses of the four planets were in a very good agreement 
with the present theoretical estimations. However, some impor- 
tant simplifications were assumed in that work. It is well known 
that planet formation, disk evolution and planetesimal migration 
occur on the same time scale. Then, a closer approach to the 
formation of the giant planets should include an approximation 
of the other two processes. In Benvenuto et al. (2009) -and also 
in Desch (2007)- the migration of planetesimals due to the gas 
drag of the disk was not taken into account. As it was shown 
by Thommes et al. (2003), Chambers (2006), and Brunini & 
Benvenuto (2008) this effect has a strong influence on the accre- 
tion time, especially for small planetesimals. Then, planetesimal 
migration should be incorporated in a more realistic model. On 
the other hand, the evolution of the gaseous component of the 
disk should be taken into account, as the gas density is an input 
parameter when determining the accretion rates and the planetes- 
imals migration velocities. Finally, the planetesimals relative ve- 
locities were considered to be in equilibrium. Nevertheless, it is 
not clear that an equilibrium condition is ever achieved, specially 
for small planetesimals (Chambers, 2006; Brunini & Benvenuto, 
2008). Instead, relative velocities out of equilibrium should be 
adopted due to the fact that the size distribution of the planetesi- 
mals extends to meter-size objects. 

These conditions were incorporated in a recent model devel- 
oped by Guilera et al. (2010). They calculated the simultaneous 
formation of Jupiter and Saturn -at their current locations- im- 
mersed in a protoplanetary disk that evolves with time. They 
found that the simultaneous formation of Jupiter and Saturn 
could be very different when compared to the isolated formation 
of each planet. For the classical solar system surface density pro- 
file (E oc cT 3 / 2 ), they found that the formation of Jupiter inhibits 
the growth of Saturn. On the other hand, for a disk with surface 
density profile I oc a~ 1/2 , Saturn turned out to form faster than 
Jupiter. Moreover, the quick formation of Saturn produced an in- 
ner solid density wave that significantly reduced the time of the 
formation of Jupiter. Finally, for a disk with surface density pro- 
file £ oc a~ l , simultaneous and isolated formation of both planets 
resulted in very similar outcomes. 

The aim of this work is to study some of the parameters that 
characterize the protoplanetary disk in order to determine the 
conditions that could favor/not favor the formation of the solar 
system giant planets, in line with the initial configuration pro- 
posed by the Nice model. The formation of the planets is cal- 
culated using the in situ approximation. However, the formation 
of the four planets are computed simultaneously, which means 
that the formation of one planet can affect the formation of an- 
other one because it perturbs the planetesimal disk. We will fo- 
cus then our analysis on the first stage of the planetary formation, 
this meaning, before the gaseous runaway starts. We aim to have 
an insight into the conditions of the disk (especially, the surface 
density profiles and the planetesimal size distribution) that lead 
to the formation of the cores of Jupiter and Saturn on a timescale 
that enable the planets to reach the gaseous runaway. Moreover, 
Uranus and Neptune should also be able to grow to their current 
mass. For this purpose, we consider a population of planetesi- 
mals with non homogeneous radii following a power law distri- 
bution, where the radii go from r™" to 100 km (the maximum 
value of the distribution is fixed for all the cases). The value of 
the minimum radius is to be determined as follows: we first cal- 
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culate the isolated formation of the four planets for different val- 
ues of y 1 ™" and determine then the optimum values they can take 
(by "optimum values" we mean those values of r^ m that the four 
planets have in common that minimize their time of formation ). 
With these values we then explore the simultaneous formation. 
This procedure is employed to analyze different profiles for the 
surface density of the protoplanetary disk. Finally, we discuss 
the scenario that best matches the basic observational constraints 
for the formation of the giant planets (the disk dissipation time 
scale and the estimated mass of the core) with the conditions of 
the Nice model. 

The reminder of this paper is organized as follows: In Section 
2 we give a brief description of our model, developed in a pre- 
vious work (Guilera et al., 2010). In Section 3 we calculate the 
formation of the giant planets of the solar system for a steep sur- 
face density profile similar to the one derived by Desch (2007). 
Then, in Section 4 we present the results for smoother surface 
density profiles. Finally, in Section 5 we discuss our results and 
present the concluding remarks. 

2. A brief description of the model 

In a previous work (Guilera et al., 2010), we developed a numer- 
ical code to compute in situ but simultaneous formation of sev- 
eral giant planets immersed in a protoplanetary disk that evolves 
with time. In the present paper, we used this code to study the 
first stage of the formation of the four giant planets of the solar 
system. 

The protoplanetary disk was characterized by gaseous and 
solid (a planetesimal disk) components. We considered a clas- 
sic power-law to describe the surface density, £ oc a~ p . In our 
model, the disk extends from 0.4 to 30 AU. We fixed the outer 
edge of the disk at 30 AU in order to be consistent with the 
outer edge of the planetesimal disk used by the Nice model. For 
the planetesimal disk, we considered a population of non-equal 
sized bodies. Planetesimals were assumed to follow a size dis- 
tribution whose radii were between a minimum value, r^ m , and 
100 km, with steps selected in order that the quotient of masses 
of consecutive sizes is a factor of two. For example, the contin- 
uous planetesimal size distribution between 1 km and 100 km 
was characterized by 21 different sizes of planetesimals, while 
the continuous planetesimal size distribution between 10 m and 
100 km was characterized by 41 different sizes of planetesimals. 
Kokubo & Ida (2000), employing N-body simulations, and more 
recently, Ormel et al. (2010), employing statistical simulations, 
showed that in the oligarchic growth regime planetesimals fol- 
low a mass distribution dnjdm oc m~ p , with p ~ 2.5. For this 
work, we adopted this distribution. Note that most of the mass is 
in small planetesimals. 

As we mentioned above, we considered the protoplanetary 
disk to evolve. In terms of the planetesimal disk, we considered 
planetesimal migration due to nebular gas drag. We adopted that 
the migration velocities are given by the prescription of Adachi 
etal. (1976), 
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Here p p is the bulk density of planetesimals (1.5 g cm 2 ), r p is 
the radius of the planetesimal, Co ~ 1, V = (vk - v gas )/vk, and 



a is the exponent of the power-law density of the gas in the disk 
mid plane (p g oc a~ a ). Then, the evolution of the planetesimals 
disk obeys a continuity equation, 



(92, 1 d I da^ 
dt ada\ dt 



F(a), 



(3) 



where S s is the surface density of solids and F(a) describes the 
sinks of disk material (accretion by the forming planets). 

Adopting a temperature profile for an optically thin disk (T oc 
a~ 1/2 ) the density distribution of gas on the mid-plane of the disk 
is given by p g oc aT p ^^ ', and for simplicity, we assumed that 
the gaseous component dissipates out following an exponential 
decay of its density with a characteristic time-scale of 6 Myr 
(Haisch et al., 2001), 

p g (a,t)=p g (a,0)e-" 6M y\ (4) 

For the growth of the core, we adopted the oligarchic regime, 
in accordance to our previous work (Guilera et al., 2010). 
Adopting the particle in a box approximation (Inaba et al., 2001) 
for the accretion of solids 



colli 



(5) 



dM c _ 2nlZ s (ap)R 2 H ^ 
dt ~ P 

and integrating over the planetesimals size distribution and the 
feeding zone, the solid accretion rate is given by 



dM c 
dt 



J dm J 

JPSD JFZ 



2mp(a,RH,ap) X 
2jiY, s {a, >ri)R 2 H 



P co u(a,m) a da. 



(6) 



where Mq is the mass of the core, Rh is the Hill radius, P is the 
orbital period, ap is the semi-major axis of the planet, and is 
a normalization function (see Brunini & Benvenuto, 2008). P co ii 
is the collision probability between the planetesimals and the 
planet, which is a function of the planet's core radius, the Hill 
radius of the planet, and the relative velocity of planetesimals 
Pcoii = Pcoii(Rc,RH,v re i) (for more details, see Guilera et al., 
2010). 

We also took into account the enhancement of the planet's 
capture cross-section due to the presence of the planet's enve- 
lope. Inaba & Ikoma (2003) found an approximate solution to 
the equations of motion, which allows for a rapid estimation of 
the radius of the planetesimal captured, r p , as function of the the 
planet's enhanced radius ^ c 



3p(R c )R c I vl + 2GMp(R c )/Rc \ 
2p p \v 2 co + 2GMp(R c )/R H )' 



(7) 



where v ro is the relative velocity of the planet and planetesimal 
when the two are far apart while M P (R C ) and p(Rc) are the total 
mass of the planet and the density of the planet's envelope con- 
tained within R c , respectively. Inaba & Ikoma (2003) propose 
replacing R c for R c in the expressions of collision probability, 

so P coU = P co u(Rc,RH,V re l). 

Planetesimal's relative velocities out of equilibrium were 
prescribed. The relative velocity, v re i, between a planetesimal 
and the protoplanet may be described by 



y,-ci = \l ~e 2 + ^' 2v k, 



(8) 



where v k is the keplerian velocity. We considered that the rel- 
ative velocity is governed by two factors, the gravitational stir- 
ring caused by the protoplanets and the gas drag damping. We 
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adopted the prescriptions developed by Ohtsuki et al. (2002) 
for the gravitational stirring and the prescriptions developed by 
Adachi et al. (1976) for the gas drag damping (see Guilera et al., 
2010, for a more detailed explanation). 

Finally, the equations governing the evolution of the gaseous 
envelope are those classic of stellar evolution theory 



dr 
dm r 
dP 
dm r 
dL r 
dm r 
dT 
dm r 



1 



Anr 2 p 
Gm r 



equation of definition of mass (9) 
equation of hydrostatic equilibrium (10) 



dS 

e„i - T— equation of energetic balance 
at 



Gm r T 



V equation of energy transport, 



(11) 
(12) 



where p is the density of the envelope, G is the universal gravita- 
tional constant, e p i is the energy release rate due to the accretion 
of planetesimals, S is the entropy per unit mass, and V = 
is the dimensionless temperature gradient, which depends on the 
type of energy transport. 

These equations were solved coupled self-consistently to the 
planetesimal's accretion rate (Eq.[6]i, employing a standard finite 
difference (Henyey) method and a detailed constitutive physics 
as described in Fortier et al. (2007; 2009), and Guilera et al. 
(2010). 

In the following sections we shall present our results for the 
simultaneous, in situ formation of the four giant planets of the 
solar system. Note that planets were assumed to be in circular 
orbits around the Sun with fixed orbital radii consistent with the 
Nice model. Then, the migration of the planets is not addressed 
in the present study. 

As we mentioned before, Thommes et al. (2008) showed that 
solar system analogues come out if gas giants do form but un- 
dergo modest migration and eccentricity growth, even for mas- 
sive protoplanetary disks. Miguel et al. (2010 a; b) did a statis- 
tical study of the formation of planetary systems and they also 
found that solar system analogues are favored in massive disks 
where there is not a large accumulation of solids in the inner re- 
gion of the disk, and only if type I migration is strongly reduced. 
Both works showed that toward low disk masses, planet forma- 
tion is too slow to produce gas giants during the disk lifetime. 

Several authors have argued that type I migration should be 
much slower than previously thought, otherwise planets get too 
close to or are even engulfed by the star (see, for instance, Alibert 
et al. (2005) for detailed calculation of the formation of Jupiter 
and Saturn including migration). Tanaka et al. (2002) derived 
the migration rate for vertically isothermal disks. However, it 
was recently found that migration rates, and even the direction 
of migration, could be very different when vertically radiative 
or adiabatic transfer for the protoplanetary disk is adopted (Kley 
et al., 2009; Paardekooper et al., 2010). In light of recent re- 
sults that suggest type I migration might not be as important as 
previously thought in more realistic disks (i.e not the idealized 
case proposed by Tanaka et al.), we find our hypothesis of in situ 
formation to be a good approximation to a more complex and 
realistic model, at least at early stages of the formation process. 

The neglect of type II migration is more problematic given 
that this is likely to be an important effect in many systems, espe- 
cially since the Nice model and its sequels calls for an early mi- 
gration and multi -resonant capture of the giant planets. However, 
calculating type II migration rates is a complicated and subtle 
procedure in multi-planet systems and it is beyond the scope of 



this paper. It is important then to make here a note of caution: 
our results should be analyzed in the context of the in situ forma- 
tion, which could be considered as a good approximation during 
the first stage of the formation of giant planets, i.e. before the 
gaseous runaway begins. 

In addition to neglecting migration, our code does not cal- 
culate the final stages of the gas accretion. Our code calcu- 
lates the gas accretion rate of the planet self-consistently using 
an adapted Henyey-type code of stellar evolution, but always 
under the assumption that the disk can supply the necessary 
amount of gas required by the planet. Then, after the runaway 
growth of the envelope sets in, this hypothesis may not be any 
more valid. However, recent works have shown that gas accre- 
tion onto a critical mass core is likely to be rapid. Machida et 
al. (2010) employing three-dimensional hydrodynamics simula- 
tions, found that the gas accretion time-scale of a giant planet 
is about 10 5 yr (two orders of magnitude less than usually esti- 
mated the time-scale to form the core). This means that the real 
bottleneck for giant planet formation is the growth of a critical 
mass core in the first place. 

For the above mentioned reasons we think our code is suit- 
able for studying the first stage of planetary formation: the 
growth of the planets up to their critical masses and how the 
simultaneous formation differs from the isolated case. 

3. Simultaneous formation of the solar system giant 
planets for a disk with a steep profile: I oc a 2 

Adopting the initial configuration proposed by the Nice model, 
Desch (2007) derived a surface density profile for an analogue 
of the minimum mass solar nebula which can be adjusted with a 
single power law, X oc cT 2168 . This profile is much steeper than 
the one derived by Hayashi (1981), where £ oc a -1 - 5 . We started 
then our investigation adopting a similar density profile as that 
proposed by Desch; in our case, and for the sake of simplicity, 
we adopted X oc cr 2 . Pollack et al. (1996) have used this profile 
to study the isolated formation of giant planets in solar system, 
but adopting a faster time-dependent regime for the growth of 
the core than the corresponding to the oligarchic growth regime. 

Following the work of Benvenuto et al. (2009), we imposed 
the value of 1 1 g cirT 2 for the initial solids surface density at the 
location of Jupiter. Then, at t = when calculations are started, 
the density profile throughout the disk is described by: 



(13) 



where r\ lce takes into account the condensation of volatiles be- 
yond the snow line, considered to be located at 2.7 AU, 



Vic, 



a > 2.7 AU, 
a < 2.7 AU. 



(14) 



In this work, we adopted a value for the gas-to-solid ratio of 
Z 1 a 65, where Z = 0.0153 (Lodders et al., 2009) is the initial 
abundance of heavy elements in the Sun, implying that the initial 
surface density of gas is given by: 



I, = 719 



(5TATj) 2gCnr 



(15) 



Rescaling to 1 AU and spreading the snow line in a region of 
about 1 AU with a smooth function as proposed by Thommes et 
al. (2003), the disk was defined as: 



83.19 + (332.75 -83.19) 



1 .la 



- tanh I - — — 
2 \ 0.5 j 
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Table 1. Isolated formation of the solar system giant planets' 
cores for a disk with surface densities of solids and gas oc gT 2 . 
Here r™ '" stands for the minimum radius of the size distribution 
of planetesimals. f r ° e s r s represents the time at which the mass of 
the envelope of the planet achieves the core's mass of the planet 
and gaseous runaway starts. M c r ° e s r s is the mass of the core at 
the cross-over time. Quantities between brackets, for the cases 
of Neptune and Uranus, correspond to the time, and the respec- 
tive core mass, at which Neptune and Uranus achieve its current 
masses (~ 17 and ~ 14.5 earth masses, respectively). When the 
cross-overtime exceeds ~ 10 Myr no values are given. 



P 

[m] 


Jupiter 

over 

[M e ] 


tcross 
l over 

[Myr] 


Saturn 

ifcross 

over 

[M e ] 


across 
l over 

[Myr] 


Neptune 

A/icross +cross 
over over 

[M e ] [Myr] 


Uranus 

Tyrcross *£ross 
over over 

[M ffi ] [Myr] 


10 






19.13 


6.96 


22.02 2.90 
(14.62) (2.31) 


18.75 6.01 
(12.59) (4.98) 


12.5 






19.22 


5.15 


21.07 3.04 
(14.57) (2.43) 


13.10 13.88 
(11.24) (12.39) 


15 






19.65 


3.91 


19.26 3.85 
(14.27) (3.17) 




20 






20.90 


2.83 


14.48 8.69 
(12.78) (7.82) 




25 


12.54 


14.98 


20.92 


2.78 






40 


16.89 


6.09 


16.30 


5.73 






50 


18.34 


4.38 


12.87 


10.43 






75 


18.80 


3.91 










100 


17.59 


4.79 










150 


15.22 


7.06 










200 


13.61 


9.43 










Table 2. Same as TableQ]but for the simultaneous formation. 


r min 
P 

[m] 


Jupiter 

Hfcross 

over 

[M ffi ] 


L over 

[Myr] 


Saturn 

Hfcross 

over 

[M e ] 


l over 

[Myr] 


Neptune 

luicross +cross 
over over 

[M ffi ] [Myr] 


Uranus 

rurcross *cross 
over l over 

[M e ] [Myr] 


10 


~ 0.5 


» 10 


17.73 


8.73 


19.81 5.16 
(14.86) (4.28) 


13.61 10.26 
(11.80) (8.24) 


50 


15.47 


8.89 


13.42 


10.44 | 


~ 1 » 10 


~ 0.5 » 10 



(lAu) 



gem 



-2 



XJa) = 21750 



(l Au) 



gem 



(16) 
(17) 



Regarding the temperature profile, we adopted the same one as 
in our previous works, 



T(a) = 280 (tau) k - 



(18) 



We did not consider the evolution in time of the temperature 
profile, which implies that the temperature of the nebula at a 
given location is a fixed external boundary condition. The initial 
mass of the disk was ~ 0.066 M Q , this value being the result 
of integrating the surface density profile from a = 0.4 AU to 
a = 30 AU. 

With the disk so defined, we first calculated the isolated for- 
mation of each of the four giant planets of the solar system. For 
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Fig. 1. Cross-over time as function of the minimum radius of the 
planetesimals size distribution for the case of the isolated forma- 
tion of each planet. Results correspond to a disk with a surface 
density profile X oc cC 1 . Open triangles (diamonds) correspond 
to the time at which Neptune (Uranus) achieves its current mass 
(~ 17 and ~ 14.5 Earth masses, respectively). 



24 



^ 20 



— 18 



16 



14 



12 



10 





Jupiter (5.5 AU) - 






Saturn (8.3 AU) - 




i Is. 


Neptune (1 1 AU) 




\ J* — •-. 




-A- 




Uranus (14 AU) 


♦ 




o 














\ * " & / 


• 











10 12.5 15 20 25 



40 50 
min 



75 100 



150 200 



[ml 



Fig. 2. Core mass at the cross-over time as function of the min- 
imum radius of the planetesimals size distribution for the case 
of the isolated formation of each planet. Results correspond to 
a disk with a surface density profile E oc a' 2 . Open triangles 
(diamonds) correspond to mass of the core at the time in which 
Neptune (Uranus) achieves its current mass (~ 17 and ~ 14.5 
Earth masses, respectively). 



each planet we run several simulations only changing the min- 
imum radius of the size distribution of planetesimals, r™". The 
aim of this procedure was to look for an interval in the plan- 
etesimals radii where the isolated formation of all the planets 
occurs in less than 10 Myr. Afterwards, and using these results 
as a guide, we looked for an optimum value of r™ to calculate 
the simultaneous formation of the four planets. 

The results for the isolated formation are summarized in 
Table [TJ We found that each planet has an independent interval 
for r""" in which the cross-over time (the time at which the mass 
of the envelope achieves the core's mass, and gaseous runaway 
starts) is reached in less than 10 Myr (Fig. [TJ. Jupiter reached 
its cross-over mass, and could be formed, in less than ~ 10 Myr 
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for r™'" in the range of ~ 30 and 200 m. Saturn, on the other 
hand, reached its cross-over mass in less than -10 Myr for r™" 
between 10 and ~ 50 m. In the case of Neptune, rj™ had to be 
between 10 and 20 m. Finally, Uranus could only be formed in 
agreement with the time restriction for r™" < 12.5 m. 

Note that we did not consider radii smaller than 10 m. This 
is due to the fact that in most cases (see next Sections) for such 
small planetesimals the time-scale of formation of the outermost 
planets are significantly lower than the time-scale corresponding 
for Jupiter and Saturn. We consider very unlikely that Uranus 
and Neptune reached their current masses before Jupiter and 
Saturn reached their cross-over masses, as in that case they could 
have continue their growth and become gas giants. 

For Desch's steep profile, we found that the shortest cross- 
over time for Jupiter is ~ 4 Myr, for Saturn is ~ 2.8 Myr, for 
Neptune is ~ 2.9 Myr and for Uranus is ~ 6 Myr. However, 
these time-scales did not occur for the same minimum radii of 
the planetesimals populating the disk. 

On the other hand, we found that for Jupiter and Saturn the 
cross-over masses were in good agreement with the current the- 
oretical estimations of their solid content (Fig. Note that we 
assumed that all the infalling planetesimals reach the core's sur- 
face without losing mass on their trajectories throughout the en- 
velope, this meaning that M c corresponds in reality to the to- 
tal heavy element's mass in the interior of the planet (core plus 
solids in the envelope). Theoretical models estimate the core 
mass of Jupiter between 0- 12 M ffl (Guillot, 2005) or 14-18 M ffi 
(Militzer et al., 2008), depending on the equation of state em- 
ployed. Regarding the core mass of Saturn, it is accepted to be 
in the range of 9 - 22 M® (Guillot, 2005). Models also predict 
10-40 M© and 20 - 30 M ffl of total heavy elements (envelope + 
core) in Jupiter and Saturn, respectively (Guillot, 2005; Guillot 
& Gautier, 2009). On the other hand, Uranus and Neptune are 
mostly rocks and ices. The planetary interior models of Podolak 
et al., (2000) place upper limits on the H/He contents of Uranus 
and Neptune, being for Uranus of 4.2 and 3.2 M ffl for Neptune. 
If only hydrogen and helium are present in the atmosphere, a 
lower limit for the gas mass of 0.5 M ffl is obtained for each planet 
(Guillot, 2005). To relate the above mentioned values to our re- 
sults of the core mass, in the following we will consider that 
the total amount of heavy elements contained in each planet cor- 
responds to the mean value of the current estimated boundaries. 
This means that we are going to compare our results for the mas s 
of the core with: 25 ± 15 M ffl for Jupiter, 25 ± 5 M e for Saturn, 
15.15 ± 1.35 M e for Neptune, and 12.15 ± 1.85 M e , for Uranus. 

Regarding the formation of Neptune and Uranus, we can see 
(Table [T]and Fig. [2]) that if we let them grow without restriction 
concerning their masses, for most cases their cross-over masses 
were higher than their current masses. However, the masses of 
the cores at the time they reached their current masses are in 
agreement with theoretical estimations. 

For the cases of Jupiter and Saturn, we found that there was 
an optimum value for rj* 1 of around 70 m for Jupiter, and be- 
tween 20 - 25 m for Saturn, for which the accretion was most 
effective. The existence of an optimum r 1 ™" for the size distribu- 
tion of the planetesimal radii is due to the fact that planetesimal's 
migration velocities are oc 1 /r p (Eq. [TJ and that the planetes- 
imal's accretion by the protoplanet is more efficient for small 
objects. Then, there is a competition between the efficiency of 
the accretion and the efficiency of the planetesimal migration. 
As we considered a planetesimal mass distribution that follows 
the power law n(m) - m~ 2 5 (where most of the mass of solids lie 
in the small objects), r™ ! " has to be chosen considering the above 
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Fig. 3. Time evolution of Jupiter's planetesimals accretion rate 
for the cases of isolated (solid line) and simultaneous formation 
(dashed line) for a disk profile E oc a~ 2 . 

mentioned compromise between the migration and the accretion 
rate. If the radii of the smaller planetesimals are very small, most 
of the mass of solids would pass through the feeding zone of the 
planet without being accreted. If, on the other hand, they are 
very large, the accretion time scale would be longer than the 
protoplanetary disk lifetime. 

The reason why the intervals for r^ m of each planet are dif- 
ferent is due to the fact that both the accretion rate of planetesi- 
mals and their migration velocities vary with the distance to the 
central star. For this steep profile, we did not find an overlap be- 
tween these intervals, so it seems that there is not a common 
size distribution for the planetesimal radii that would allow us to 
form simultaneously the four planets in less than 10 Myr. 

However, in spite of the previous results, we checked if ef- 
fectively there was not an interval where the simultaneous for- 
mation could occur. Taking into account the results for the calcu- 
lations of isolated formation, we run some simulations where we 
considered the simultaneous formation of the four giant planets. 

Results are listed in Table [2] We first chose r™"' = 50 m 
for the size distribution of planetesimals. In this case the for- 
mation times for Jupiter and Saturn were of ~ 9 Myr and 
~ 10.5 Myr, respectively. However, in this time interval the 
embryos of Neptune and Uranus achieved only a core mass of 
M c ~ 1 M ffi and M c ~ 0.5 M ffi respectively, both with a negligible 
envelope (the simulation was stopped at ~ 1 1 Myr, when Saturn 
finished its formation). We note that the presence of Saturn de- 
layed the formation time of Jupiter by a factor ~ 2 when we 
compare this results to the isolated ones. As we can see in Fig. [3] 
the presence of Saturn significantly decreased the solid accre- 
tion rate of Jupiter. The reduction in the solid accretion rate is 
due to the fact that planetesimals that were accreted by Jupiter 
in the isolated case (coming from the outer region of the Solar 
System), were accreted first by Saturn when the formation of 
both planets was calculated simultaneously. For this steep pro- 
file, the presence of Saturn acts as a sink of planetesimals, signif- 
icantly decreasing the surface density of smaller planetesimals 
in the neighborhood of Jupiter when we compare it with the iso- 
lated formation of Jupiter (see Fig. |4}. 

When we considered smaller planetesimals (r™ m = 10 m), 
Neptune turned out to be the first one to achieve its cross-over 
mass, followed by Saturn and Uranus, while Jupiter was not able 
to reach its cross-over mass in less than 10 Myr (the Simula- 
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Fig. 4. Time evolution of the surface density of solids for small 
planetesimals in a disk with power index p = 2. Thin lines corre- 
spond to the isolated formation of Jupiter, thick lines correspond 
to the simultaneous formation of the giant planets. We can see 
how the presence of Saturn significantly reduces the surface den- 
sity of solids at Jupiter's location by the accretion of small plan- 
etesimals. The presence of Saturn delayed the formation time of 
Jupiter by a factor ~ 2 (color figure only available in the elec- 
tronic version). 



tion was stopped at ~ 15 Myr and Jupiter's core mass was only 
~ 0.5 M ffl with almost no gas bound). However, in the simulta- 
neous case, the formation of all the planets took longer than in 
the isolated case. We consider that a case like this is very un- 
likely because, while Jupiter could not grow up to its cross-over 
mass, Neptune reached its cross-over mass before Saturn, which 
means that it could have started the gaseous runaway (and be- 
come a gas giant) before Saturn completed its formation. While 
it is widely accepted that Uranus and Neptune did not achieve 
the gaseous runaway to become gas giants, we remark that in our 
simulations we did not halt the accretion onto them when they 
reach their current masses and we allowed Neptune and Uranus 
to grow up to the cross-over mass, as the perturbations that a 
planet produces on the planetesimal disk, and how it affects the 
formation of the others, depends on its mass. This is in line with 
our working hypothesis of studying the formation of the plan- 
ets up to their cross-over masses, a process that is well described 
with our code. In this work we do not attempt to explain the final 
masses of the planets. 



Then, for a nebula like the one considered here we would 
need that r™" decreases with the distance from the central star 
in order to be able to form the four planets simultaneously. 
Moreover, an accurate treatment demands that the simultaneous 
formation should be recalculated with a non homogeneous dis- 
tribution. However, in principle, there's no physical reason to 
support such a decrease in the minimum radii of the planetesi- 
mals as we move away from the Sun. In fact, a decrease in the 
planetesimals' radii would be in contradiction with the results 
found by Chambers (2010), where planetesimal size should in- 
crease with distance. 



4. Simultaneous formation of the solar system giant 
planets for shallow disk density profiles 

In the previous Section we showed that with a steep profile, as 
the one derived by Desch (2007), for a minimum mass solar neb- 
ula compatible with the Nice model, there is not a minimum ra- 
dius for the population of the planetesimals that favors the simul- 
taneous formation of the four giant planets of the solar system. 
In this section we will analyze other density profiles, all of them 
power laws where we change the value of the index p. 

4.1. The case of Hoc a~ 3/2 

Employing the standard surface density profile S °c er 3 ^ 2 we 
repeated the calculations of Sec. [3] The same procedure as in 
the previous section was employed: we calculated first the iso- 
lated formation in order to look for candidate values for r™'" and, 
afterwards, with this values we calculated the simultaneous for- 
mation of the four planets. 

Imposing again the value of 1 1 g cm 2 for the solid surface 
density at the location of Jupiter, using the same solid/gas ratio, 
and rescaling to 1 AU, the disk is defined by 



!ja) = { 35.5 + (142-35.5) 

- 3 / 2 



1 a -2.7 
- tanh 

2 \ 0.5 
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+ - 
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(l Au) 



gem 

-3/2 



gem 



(19) 
(20) 



According to this definition, the initial mass of the disk was 
~ 0.063 M . 

The results obtained for the isolated formation of each planet 
are summarized in Table [3] We note that these results are quali- 
tatively similar to those obtained for the steeper profile 2 oc cC 1 . 
However, in this case, there are some values for r™" (between 
20 - 40 m) that may allow the simultaneous formation of all the 
planets with a common size distribution of planetesimals. 

For the case of the isolated formation of each planet, Jupiter 
reached its cross-over mass in less than 10 Myr for 20 m < r™" <; 

< 160 m, Neptune for 10 m < 



400 m, Saturn for 10 m < r% m 

rf" < 80 m, and Uranus for 10 m < r™'" 
this protoplanetary nebula, there should be an optimum radius 
r mm f or w hj c h the four planets achieved their cross-over masses 
in less than 10 Myr (Pig. [5j, this radius lying in the interval [20, 

401 m. Note that for a planetesimal's size distribution between 
r mm _ 20 m and r™ flv = 100 km, the first planet to reach its cross- 
over mass was Neptune, followed by Saturn, Uranus and Jupiter 
(see Table [3}. However, for a size distribution of planetesimals 
between r™" = 35 m and r™ flv =100 km, and between tf n =40m 
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Table 3. Same as Table Q]but for a disk with surface density of 
solids and gas oc at 3 / 2 . 



P 

[m] 


Jupiter 

jtyjcross jcross 
over over 

[Me] [Myr] 


Saturn 

xjcross across 
over over 

[M e ] [Myr] 


Neptune 

Ateross jcross 
over over 

[M e ] [Myr] 


Uranus 

jurcross +cross 
over over 

[M e ] [Myr] 


10 




23.49 4.14 


28.28 1.10 
(15.30) (0.75) 


32.20 0.61 
(13.67) (0.35) 


20 


19.88 9.85 


26.93 1.33 


30.33 0.87 
(15.58) (0.61) 


26.15 1.91 
(13.37) (1.45) 


35 


20.92 4.58 


28.38 1.22 


24.02 2.52 
(15.05) (2.08) 


17.57 7.62 
(12.45) (6.50) 


40 


21.65 3.76 


27.53 1.44 


22.33 3.28 
(14.85) (2.77) 


15.19 10.92 
(11.90) (9.56) 


80 


23.67 2.61 


21.08 3.90 


13.11 13.20 
(12.3) (12.6) 





160 | 21.03 3.70 | 14.59 10.73 | — — | — — 
400 I 15.98 7.81 I - — I — — I — — 



Table 4. Same as Table [3] but for the simultaneous formation of 
the four planets. 
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Fig. 5. Same as Fig.[TJbut considering a disk with a surface den- 
sity profile S oc a~ 3 / 2 . The gray zone corresponds to the interval 
of values of the minimum radius of the size distribution of plan- 
etesimals where the cross-over times corresponding to the four 
planets are less than ~ 10 Myr. Squares, circles, triangles and 
diamonds correspond to Jupiter, Saturn, Neptune and Uranus, 
respectively. 



V 

[m] 


Jupiter 

A/rcross jcross 
over over 

[M e ] [Myr] 


Saturn 

Jurcross ftrross 
lvl over L over 

[Me] [Myr] 


Neptune 

A/tcross *cross 
over over 

[M ffi ] [Myr] 


Uranus 

A/rcross *cross 
over over 

[M e ] [Myr] 


20 


35.74 0.99 


38.30 0.55 


33.65 0.60 
(16.20) (0.42) 


27.83 * 1.6 * 
(13.92) (1.1) 


35 


31.89 2.25 


29.12 1.83 


23.43 3.88 
(15.52) (3.01) 


16.74 8.68 
(12.69) (6.80) 


40 


30.08 2.96 


27.28 2.55 


20.72 5.27 
(15.12) (4.19) 


14.70 12.23 
(12.18) (10.02) 



* For this case, the numerical models corresponding to Uranus did 
not converge (Uranus did not achieve its cross-over mass) and the 
simulation was halted. At time t = 1.6 Myr., Uranus reached a core 
mass of 27.83 M ffi with a corresponding envelope mass of 8.75 M e . 



and r"p lx = 100 km, the first planet to reach the cross-over mass 
was Saturn, followed by Neptune, Jupiter and Uranus. 

We calculated then the simultaneous formation of the four 
planets for r™" = 20, 35 and 40 m. Results are summarized in 
Tableland Figure [8] Note that for i 3 f n = 20 m the formation 
of the system changed drastically. The rapid formation of the 
outermost planets (especially the rapid formation of Saturn) sig- 
nificantly favored the formation of Jupiter. Similar results were 
found in our previous work (Guilera et al., 2010), in the cases 
where Saturn was formed before Jupiter. This rapid formation 
of the outermost planets induced an inner density wave respon- 
sible for the reduction in the time of formation of Jupiter and 
Saturn. The density wave (see Fig.|6]and Fig. [TJ is a product of 
the excitation of the planets over the planetesimals. The gravita- 
tional stirring of the planets increases planetesimal's eccentricity 
and inclination, and force the inner migration of planetesimals 
(see Eq. [TJ, enhancing the solids surface density and, this way, 
speeding up the formation of the cores of the planets. We note 
that we only considered the simultaneous formation of four em- 
bryos. However, oligarchic growth predicts the presence of sev- 



eral embryos separated by ~ 10 mutual Hill radii. The presence 
of several embryos may damp the solid density wave, but on the 
other hand, may also favor the formation of massive cores by 
their mutual collisions, which may lead to their fusion. 

For r"™ = 35 and 40 m, the formation of the system was 
quantitatively different. Although the formation time of Saturn 
was increased by the presence of the other embryos, it could 
still be considered fast. The rapid formation of Saturn signifi- 
cantly decreased the formation time of Jupiter. As we can see in 
Fig - IS] and Fig. [7] the quick formation of Saturn forced the mi- 
gration of planetesimals and augmented the surface density of 
solids in Jupiter's neighborhood. After 0.5 Myr. the surface den- 
sity of solids in Jupiter's neighborhood was increased for small 
planetesimals (Pig. |6j due to the presence of Saturn. For bigger 
planetesimals, the shape of the density wave is more evident, and 
in some cases the value of the surface density of solids become 
greater than the initial ones (Pig. |7J. 

In these two simulations, the cross-over times and cross-over 
masses of Jupiter and Saturn were in good agreement with obser- 
vation and theoretical estimations of protoplanetary disks life- 
times and current solids content in the interior of the planets re- 
spectively (Fig. [8). Moreover, the time-scale at which both plan- 
ets achieved the gaseous runaway were very similar. Regarding 
Uranus and Neptune, their cross-over masses were greater than 
their current total masses. However, for both planets, the core 
masses at the time they reached their current masses were in 
agreement with theoretical estimations. We remark that the for- 
mation time-scale of the gas giants was shorter than the forma- 
tion time-scale of the ice giant. 



4.2. The case of Hoc a 1 

In this section, we repeat the previous methodology but now em- 
ploying a surface density profile for the disk of the form S oc cC l . 
Then, the disk is defined by 



Y s (a) 



15.125 + (60.5 - 15.125) 
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Fig. 6. Time evolution of the surface density of solids for small 
planetesimals at Jupiter-Saturn region, for a disk with a power 
index p = 3/2, considering the case of the simultaneous for- 
mation and employing a distribution of planetesimals between 
r™'" = 35 and r™ ax = 100 km. After 0.5 Myr. the surface density 
of solids at Jupiter's neighborhood is significantly augmented by 
the incoming flow of planetesimals due to the gravitational per- 
turbations of the outer planets (color figure only available in the 
electronic version). 



\1 AU7 



gem 



S„(a) = 3954.25 



(l Au) 



g cm 



-2 



(21) 
(22) 



The initial mass of the disk is ~ 0.082 M . 

In Table [5] and Fig. [9] we show the results corresponding to 
the isolated formation of each planet. Jupiter reached its cross- 
over mass in less than 10 Myr for 10 m < r™" < 800 m, Saturn 
for 10 m < r'p'" < 400 m, while Neptune and Uranus reached 
their cross-over mass in less than 10 Myr for 10 m < r™ m < 
200 m. Then, for r"™ 1 between 10 m and 200 m the isolated for- 
mation of the four planets occurred in less than 10 Myr. 

According to these results we calculated the simultaneous 
formation of the giant planets for a size distribution of planetesi- 
mals between r™'" and 100 km, where we adopted for r™" several 
discrete values between 10 m and 200 m. 

For the size distributions of planetesimals between 10 m- 
100 km and 50 m - 100 km, the rapid formation of the outer 
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Fig. 7. Same as Fig. |6]but for planetesimals with radii larger than 
1 km. In this case the shape of the density wave is well defined. 
The density wave increases the surface density of solids to val- 
ues that are even larger than the initial ones (color figure only 
available in the electronic version). 



planets significantly augmented the planetesimals accretion rate 
of Jupiter and Saturn. Planetesimals accretion rate became so 
high (10~ 3 - 10~ 2 M e yr~') that models corresponding to Jupiter 
and Saturn did not converge. We note that the mass of the cores 
corresponding to Jupiter and Saturn, before models did not con- 
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Fig. 8. Cross-over time (top) and cross-over mass (bottom) as 
function of the minimum radius of the size distribution of plan- 
etesimals for the case of the simultaneous formation. Here the 
disk surface density is characterized by S oc cT 3 / 2 . For r™" = 35 
and 40 m the cross-over times and cross-over masses of Jupiter 
and Saturn are in good agreement with the observational and the- 
oretical estimations. Regarding Neptune and Uranus their cross- 
over masses are higher than their current total masses. However, 
for both planets the mass of the core at the time they reach their 
current masses (open triangles and diamonds) are in agreement 
with theoretical estimations. 



verge, were larger than the total heavy element corresponding to 
Jupiter (~ 40 M e ) and Saturn (~ 30 M ffi ), respectively. 

On the other hand, for rf n = 100, 150 and 200 m, the re- 
sults for the simultaneous formation were not very different from 
the isolated ones (see Table [SJ. We remark that for these cases 
our results for both the cross-over times and cross-over masses 
of Jupiter and Saturn were in nice agreement with the observa- 
tions and theoretical estimations, respectively (Fig. [TOV Again, 
the cross-over masses of Uranus and Neptune were larger than 
their current masses, but the core masses at the time they reached 
their current masses were in agreement with theoretical estima- 
tions. We note again the important fact that the time-scale of 
formation of the gas giants was shorter than the corresponding 
one to the ice giants. 



Table 5. Isolated formation of the solar system giant planets' 
cores for a disk with a surface density of solids and gas oc a . 



jjnin 
P 

[m] 


Jupiter 

A/icross fcross 
over over 

[Me] [Myr] 


Saturn 

jurcross across 
over over 

[M ffi ] [Myr] 


Neptune 

njcross across 
over over 

[M e ] [Myr] 


Uranus 

A/icross fCTOSS 
over over 

[M e ] [Myr] 


10 


24.66 10.49 


25.38 3.00 


30.71 0.59 
(15.48) (0.20) 


42.56 0.09 
(14.20) (0.05) 


50 


26.32 2.31 


34.41 0.68 


32.81 0.85 
(15.75) (0.63) 


30.64 1.19 
(13.61) (0.90) 


100 


27.46 2.18 


28.39 1.84 


26.20 2.55 
(15.31) (2.19) 


23.55 3.80 
(13.17) (3.21) 


150 


25.99 2.47 


25.16 2.82 


22.51 4.36 
(14.87) (3.85) 


19.52 7.09 
(12.72) (6.16) 


200 


25.00 2.82 


22.94 3.85 


19.85 6.45 
(14.40) (5.81) 


16.48 11.49 
(12.20) (10.25) 



400 | 21.06 4.61 | 17.28 9.11 | — — | — — 
800 I 16.16 9.78 I — - I — — I — - 



Table 6. Same as Table |5]but for the simultaneous formation of 
the giant planets. 



yinin 
P 

[m] 


Jupiter 

M c tf 

[M e ] [Myr] 


Saturn 

M c tf 
[M ffi ] [Myr] 


Neptune 

M c tf 

[M e ] [Myr] 


Uranus 

M c tf 

[M e ] [Myr] 


100 


32.03 2.14 


28.53 1.98 


26.97 2.81 
(15.61) (2.27) 


24.25 4.18 
(13.50) (3.13) 


150 


28.87 2.73 


23.85 2.80 


22.01 4.97 
(15.23) (4.11) 


19.08 7.15 
(12.90) (5.90) 


200 


26.54 2.96 


19.17 3.82 


18.83 7.68 
(14.65) (6.43) 


16.15 11.26 
(12.40) (9.55) 




150 200 



800 



[m] 



Fig. 9. Same as Fig.|5]but for a disk with a surface density profile 
2 oc a , 



4.3. The case of 2 oc a 



-1/2 



Finally, we repeated the calculations employing a very smooth 
surface density profile, proportional to a In this case, the 
resulting expressions for £ are: 
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Fig. 10. Cross-over time (top) and cross-over mass (bottom) as 
function of the minimum radius of the size distribution of plan- 
etesimals for the case of the simultaneous formation. Here the 
disk surface density is characterized by X oc aT x . For all cases the 
cross-overtimes and cross-over masses of Jupiter and Saturn are 
in good agreement with the observational and theoretical estima- 
tions. Regarding Neptune and Uranus their cross-over masses 
are higher than their current total masses. However, for both 
planets the mass of the core at the time they reach their current 
masses (open triangles and diamonds) are in agreement with the- 
oretical estimations. 



6.45 + (25.8 - 6.45) 

-!/2 



I—) 

\ 1 AU/ 

S g (a) = 1686.25 (— ) g 



1 (a-2.1\ 1 
2 tanh (-05 



gem 

-1/2 



-2 



(23) 
(24) 



The initial mass of the disk in this case was ~ 0.13 M . 

In Table [7] we show the results corresponding to the isolated 
formation of each planet. We see that each planet reached their 
cross-over mass in less than ~ 10 Myr with a size distribution 
of planetesimals ranging from t"™ to 100 km, where r""" could 
take values between 10 m and 1 km (Fig.fTTI). However, for this 
density profile, there were cases where the final solid masses in 
the planets interior were higher than the theoretical estimations. 



Table 7. Isolated formation of the giant planets for a disk with a 
surface density of solids and gas oc cr 1 / 2 . 



P 

[m] 


Jupiter 
M c t f 
[Me] [Myr] 


Saturn 

M c t f 
[M ffi ] [Myr] 


Neptune 

M c t f 
[M e ] [Myr] 


Uranus 
M c t f 
[M e ] [Myr] 


10 


25.27 7.62 


26.00 2.33 


34.55 0.09 
(16.50) (0.03) 


59.85 0.03 
(14.30) (0.02) 


50 


28.65 2.14 


42.00 0.10 


54.15 0.06 
(16.52) (0.04) 


55.50 0.06 
(14.20) (0.04) 


100 


31.31 2.02 


35.00 0.77 


36.73 0.44 
(15.99) (0.30) 


38.21 0.30 
(13.95) (0.17) 


500 


23.86 4.19 


23.65 4.42 


23.32 4.53 
(14.95) (4.03) 


23.38 4.00 
(13.14) (3.37) 


1000 


18.64 8.73 


17.88 10.14 


17.57 10.65 
(13.81) (9.86) 


17.83 9.72 
(12.42) (8.60) 


Table 8. Same as Table|7]but for the simultaneous formation. 


P 

[m] 


Jupiter 

M c t f 
[M e ] [Myr] 


Saturn 

M c t f 
[M e ] [Myr] 


Neptune 
M c t f 
[M e ] [Myr] 


Uranus 
M c t f 
[Me] [Myr] 


500 


53.37 2.42 


35.54 2.37 


22.62 2.92 
(15.67) (2.29) 


23.20 4.50 
(13.04) (3.78) 


1000 


47.90 4.80 


31.05 4.73 


16.85 6.07 
(13.87) (5.06) 


16.21 11.50 
(12.24) (10.09) 



The results of the simultaneous formation for this profile are 
summarized in Table [8] 

For r^ m between 10 and 100 m, the fast formation of 
Neptune and Uranus significantly augmented the planetesimals 
accretion rate of Jupiter and Saturn. Planetesimals accretion rate 
became so high (10~ 3 - 10~ 2 M ffi yr _1 ) that models correspond- 
ing to Jupiter and Saturn did not converge. We note again that the 
mass of the cores corresponding to Jupiter and Saturn, before 
models did not converge, were larger than the total heavy ele- 
ment corresponding to Jupiter (~ 40 M e ) and Saturn (~ 30 M e ), 
respectively. For r™'" = 500 m and r 11 ™ = 1 km the presence 
of Uranus reduced the formation time of the inner planets by 
a factor of 2, but increased significantly its final solid masses. 
While the cross-over time of Jupiter and Saturn are in nice agree- 
ment with theoretical estimations, the cross-over mass corre- 
sponding to Jupiter and Saturn are larger than the total heavy 
elements content corresponding to Jupiter (~ 40 M ffi ) and Saturn 
(~30 M e ). 



5. Discussion and conclusions 

In this paper we studied the in situ, simultaneous formation of 
the Solar System giant planets in the framework of the core ac- 
cretion model, and according to the conditions imposed by the 
Nice model. The goal of this work was to delimit some parame- 
ters that describe the protoplanetary disk (especially the surface 
density profile and planetesimal's population size distribution) 
in relation to the likelihood of the formation of the four giant 
planet's cores, in less than ~ 10 Myr and with a content of heavy 
elements in good agreement with current estimations. 

The surface density profile of the disk was assumed as a 
power law, X oc aT 1 ' ', where the power index p is considered 
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Fig. 11. Same as Fig. [5] and Fig. [9] but for a disk with a surface 
density profile E oc a ^ 2 . 



as a free parameter, and several values of p were analyzed. We 
considered that accreted planetesimals follow a power law mass 
distribution of the form n(m) oc m~ , where most of the mass 
is in smaller objects. We assumed planetesimals to be spheres 
of constant density, which in turn implies a power law distribu- 
tion for planetesimal's radii. We discretized the continuous plan- 
etesimal population, considering that the minimum radius, r™'", 
is also a free parameter and the maximum radius is fixed in 

jmax = 100 km. 

We note that, although we calculate the evolution of the plan- 
etesimal disk due to planet's accretion and planetesimal migra- 
tion, we do not calculate the planetesimal's evolution due to 
aggregation and fragmentation. Inaba et al. (2003) found that 
large amounts of mass could be lost due to inward drift of 
small collision fragments. So, fragmentation seems to play an 
important role and it should be considered in a more accurate 
model. Although we can not predict how fragmentation could 
affect the results of our model, we want to point out some dif- 
ferences between the working hypothesis of Inaba et al. (2003) 
and the one presented here. Inaba et al. (2003) started with an 
homogeneous population of planetesimals of radii ~ 10 km. 
So, smaller planetesimals appear due to fragmentation of big- 
ger ones. This implies that most of the mass relies in bigger 
planetesimals (Dohnanyi, 1969; Wetherill & Stewart, 1993). On 
the other hand, we started our simulations already in the oli- 
garchic growth regime, with an embryo with the mass of the 
Moon (~ 0.01 M e ) immersed in a swarm of planetesimals that 
follow a mass distribution. Ormel et al. (2010) employing statis- 
tical simulations that include several physical processes such as 
dynamical friction, viscous stirring, gas drag, and fragmentation, 
found that starting with an homogeneous populations of plan- 
etesimals of radii ro, the transition between the runaway growth 
and oligarchic growth is characterized by a power-law size dis- 
tribution of mass index q ~ 2.5, where most of the mass lie 
in small planetesimals. So, our initial conditions are consistent 
with the oligarchic growth regime using ro < 1 km. Another im- 
portant difference is the moment at which the gas drag of the 
envelope become efficient. Inaba et al. (2003) found that when 
the mass of the embryo reach ~ the mass of Mars it is capable 
to acquire an envelope. However, we found that when the mass 
of the embryos is ~ 0.1 M ffi it already has an envelope and the 
envelope gas drag becomes efficient for smaller planetesimals 
(tens and hundreds meter-size planetesimals). Moreover, as we 



can see in Fig. [T2l when the mass of the embryo is ~ 0.6 M e 
the ratio between the enhanced radius and the core radius is a 
factor ~ 7. The enhanced of the cross-section capture due to the 
presence of the envelope and the moment when it become effi- 
cient are very important. As an example, we calculated the iso- 
lated formation of Jupiter for a disk with density profiles oc a -1 , 
and with a planetesimal size distribution between 100 m and 
100 km, but now not taking into account the enhancement of 
the capture cross-section due to the envelope gas drag. We com- 
pared this simulation to the one with the same initial conditions 
but considering the enhancement of the capture cross-section. In 
Fig. [13] we show the comparison between the two simulations. 
For the case where the enhanced capture cross-section is consid- 
ered, Jupiter reached its cross-over mass at 2.18 Myr. (the value 
of the cross-over mass is 27.46 M ffl , see Table 0, while if the 
enhanced capture cross-section is not considered, Jupiter did not 
reach its cross-over mass (the simulation was halted at 15 Myr. 
and Jupiter only reached a core of ~ 2.5 M e with a negligible 
envelope). We see that when the core -for the case where the 
enhancement of the capture cross-section is not incorporated- 
reached the mass of Mars, the corresponding core for the case 
where the enhanced capture cross-section is considered, is able 
to grow up to 10 M e in the same elapsed time. It is clear that the 
enhancement of the cross-section of the core due to the presence 
of the envelope plays a fundamental role in the formation of a gi- 
ant planet. We also can see how the planetesimal accretion rate 
suffers a significant drop when the enhancement of the cross- 
section is not taken into account. Finally, we can see that the dif- 
ference between both simulation begins when the envelope gas 
drag become effective for small planetesimals (at ~ 0.01 Myr, 
which correspond to a core mass of ~ 0.1 M e ). So, the moment 
(or the mass of the embryo) at which the envelope gas drag be- 
comes effective for smaller planetesimals is important due to the 
fact that we are considering that most of the mass lies in small 
objects. The difference in the ratio between the enhanced radius 
and the core radius, respect to Inaba et al. (2003), may corre- 
spond to the fact that we are employing a more accurate model 
to calculate the growth and evolution of the planet's envelope. 
Therefore, we conclude that fragmentation could be an impor- 
tant effect, and we plan to include it in future improvements of 
our model, but we think that a definitive evaluation of its rel- 
evance should await calculations including a realistic model of 
the process, combined with a detailed prescription of the effect 
of the planet atmosphere on the fragments. 

We started studying a Desch-like profile for the disk surface 
density, S oc aT 1 (see Desch, 2007). We first computed the iso- 
lated formation of the four giant planets looking for a common 
interval for the minimum radius of the planetesimal population 
that could allow the simultaneous formation of all the planets. 
However, we did not find overlapping values for the four planets. 
Despite that, we anyway calculated some simulations in order to 
check if the simultaneous formation was possible (see Table 13, 
and confirmed that, according to our model, adopting Desch pro- 
file makes the formation of the giant planets very unlikely, as the 
formation time-scale largely exceeds the 10 Myr barrier. 

We proceeded then with less steep profiles. We analyzed 
other three cases: S oc a -1 - 5 , "L oc aT x and Z oc cr 5 , the first one 
corresponding to the standard minimum mass solar nebula and 
the other two following results of Andrews et al. (2009; 2010), 
which are based on observational data of circumstellar disks. In 
all these cases, when first calculating the isolated formation of 
Jupiter, Saturn, Neptune and Uranus, we were able to find sev- 
eral common values for r™" that could lead to a satisfactory re- 
sult when calculating the simultaneous formation (see Figs. [5] 
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M core [M e ] 

Fig. 12. The ratio between the enhanced radius and the core ra- 
dius as function of the core mass for the isolated formation of 
Jupiter, for a disk with a density profile oc a -1 , and with a plan- 
etesimal size distribution between 100 m and 100 km. The en- 
velope gas drag for smaller planetesimals becomes effective at 
very small core masses. 

[9] fTTT i. Indeed, for these values we run the corresponding sim- 
ulations for the simultaneous formation. For these three density 
profiles, the simultaneous formation could be, in principle, con- 
sidered as possible if the most abundant planetesimals are those 
whose radii range from several tens to several hundreds of me- 
ters depending on the case, but not larger than that. From our 
results it is clear that the smoother the profile density of the disk, 
the more efficient the formation process. Moreover, considering 
lower values for p allowed us to increase the minimum radius 
of the accreted planetesimals. We note that decreasing p is di- 
rectly related to considering more massive nebulae, which in turn 
means that the mass in the protoplanets feeding zone is larger, 
favoring this way the accretion process. Furthermore, planetes- 
imal migration is less efficient because gas drag is weaker, so 
planetesimals remain longer in the feeding zone of the proto- 
planets. Beside this, for a solid surface density profile E s oc a~ p 
where p < 1 the mass of solids grows outwards the planetesimal 
disk, and the incoming mass flux in the feeding zone of the pro- 
toplanets is greater than the out-coming one, which also favors 
the accretion. 

On the other hand, Desch's profile is compatible with a de- 
cretion disk. This kind of profiles turned out to be almost sta- 
tionary for about ten million years. Crida (2009) demonstrated 
that in a disk with this density profile, the four giant planets of 
the solar system would have been unable to survive. Particularly, 
Jupiter would have become a hot giant planet. Morbidelli & 
Crida (2007) and Crida (2009) also showed that Jupiter and 
Saturn could avoid migration if a less massive nebula, like 
Hayashi's nebula is assumed. While the initial value of the gas 
density used in our simulations, at the position of Jupiter, is 
about 5 times higher than that of Hayashi's minimum mass solar 
nebula, this value is decreased by an exponential dissipation fac- 
tor during the formation of the planets. Then, by the time Jupiter 
and Saturn finished their formation, the disk conditions corre- 
spond to a much less massive disk, especially in the cases where 
the formation times of Jupiter and Saturn are similar to the dis- 
sipation time-scale of the gas disk. This should help to prevent 
Jupiter and Saturn migration in initially massive protoplanetary 
disks. 



We found several values of r™" for each profile considered 
here (except p = 2), for which the simultaneous formation of the 
planets could be achieved in less than 10 Myr and the amount 
of heavy elements in the interior of the planets were in good 
agreement with theoretical estimations of current abundances. In 
most of the cases, the formation time-scale of the gas giants was 
shorter than the formation time-scale of the ice giants. However, 
when very small planetesimals were the most abundant bodies in 
the disk, the icy giants formed before the gas giants. We consider 
these cases very unlikely, at least under the hypothesis of our 
model, as they would need extra explanations to justify why they 
did not started the gaseous runaway before Jupiter and Saturn if, 
in principle, there was still plenty of material that could be ac- 
creted. This specially applies for the cases where Uranus and/or 
Neptune were the first to be formed. 

There is another interesting result we would like to high- 
light. When we consider S oc cT x and r""" = 100 m (Sec. 14.21 ) 
the cross-overtimes of Jupiter and Saturn and the time at which 
Neptune and Uranus reached their current masses are of the same 
time-scale; around 2 -3 Myr (see Tableland FigfTUt. This time- 
scale is shorter than the dissipation time-scale of the disk. The 
fact that the four planets formed on the same time-scale, which 
is also comfortably short, is not a minor point, but on the con- 
trary it presents a scenario that should be explored in more de- 
tail in the future. If the gaseous component of the disk could be 
rapidly evaporated by some external mechanism [for instance, 
if the solar system is strongly irradiated by an OB stellar as- 
sociation (see, for example, the work by Clarke (2007))], our 
results can provide a pathway to find suitable parameters to de- 
scribe the protoplanetary disk. Then, we note that there are cases 
where calculations of simultaneous formation with an appropri- 
ate planetesimal size distribution, and together with a smooth 
surface density profile, provide suitable conditions for the timely 
formation of the external planets, the four formation processes 
occurring on the same time-scale. 

Finally, we would like to remark the need that the size distri- 
bution of planetesimals extends to objects with radii < 1 km to 
form the four planets on a time-scale compatible with the obser- 
vational estimates. It was already pointed out in Goldreich et al. 
(2004) that the presence of a large amount of small planetesimals 
would help to speed up planetary formation. If the planetesimal 
disk is dominated by big bodies the time-scale to form solid em- 
bryos able to bind a significant envelope and start the gaseous 
runaway accretion would be too long to complete the formation 
in less than 10 Myr. 

However, this result contradicts recent models for planetes- 
imal formation. Johansen et al. (2007), Cuzzi et al. (2008) and 
Youdin (2011) describe different mechanisms for planetesimal 
formation incorporating turbulence as a way to aid planetesi- 
mals to grow, showing that, if planetesimals formed this way, 
they are likely to be large (~ 100 km, or even bigger, grow- 
ing directly from sub-m particles). However, these models do 
not predict a specific planetesimal size distribution and the na- 
ture of the turbulence is controversial. Nelson & Gressel (2010) 
showed that fully developed magneto-hydrodynamic turbulence 
in protoplanetary discs would have a destructive effect on em- 
bedded planetesimals. Arguments in favor of the hypothesis of 
primitive initial large planetesimals are found in Morbidelli et 
al. (2009), where through models of the collisional evolution of 
the asteroid belt they conclude that the primitive asteroids should 
have been big. However, their model does not take into account 
a very important effect: the primitive intense bombardment of 
the asteroid belt by outer planetary region comets (Gil Hutton 
& Brunini, 1999) that should change this conclusion. Moreover, 
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Weidenschilling (2010) proposed an alternative scenario to the 
work of Morbidelli et al. (2009). He found that coagulation from 
small planetesimals (~ 100 m of diameter) represents better the 
size distribution of the asteroid belt. Such small planetesimals 
could be formed by coagulation in collision driven processes by 
size-dependent drift due to nebular gas drag (Weidenschilling, 
1997). For these reasons, we consider that the planetesimal size 
issue is far from being fully understood, and a distribution where 
most of the solid mass accreted by a protoplanet comes from 
small planetesimals cannot be ruled out. 

However, one possible way to relax this assumption is to 
consider giant impacts and mergers during the planetary for- 
mation (Li et al., 2010). If a forming giant planet is impacted 
by a Mars- to an Earth-mass embryo, the impactor most likely 
reaches the core possibly remixing the core into the envelope 
(Li et al., 2010). Broeg & Benz (2011) have studied the effect 
of this impact on the gas accretion rate: initially, most of the en- 
velope can be ejected, but afterwards gas is reaccreted very fast 
and the overall gas accretion rate turns out to be larger than in 
the standard, continuous-planetesimal-accretion scenario. As a 
consequence, the planet gains large amounts of mass in a short 
time-scale and the growth-time-scale is reduced. 

Another path to accelerate the formation of the gaseous plan- 
ets is to consider the fusion of embryos as a mechanism to ob- 
tain massive cores, especially in the early stages of the disk 
evolution. Oligarchic growth predicts the simultaneous forma- 
tion of many embryos on orbits separated by about 10 Hill radii 
from each other. Some works have been done to study the fu- 
sion of growing protoplanets but only taking into account their 
solid cores, neglecting the presence of the envelope (Chambers, 
2006; Brunini & Benvenuto, 2008). To perform a full-detailed 
simulation in the context of the present study we would need 
to compute the merger of embryos self-consistently, taking into 
account the growing gaseous envelope. This is a very complex 
phenomenon, beyond the scope of the present paper, but that will 
be investigated in a future work. 
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Fig. 13. Top and middle: Time evolution of the core mass and the 
planetesimal accretion rate for the isolated formation of Jupiter 
employing a disk with density profiles oc a -1 , and a planetesimal 
size distribution between 100 m and 100 km, with (solid line) 
and without (dashed line) considering the enhancement of the 
capture cross-section due to the envelope gas drag. Bottom: The 
ratio between the enhanced radius and the core radius as function 
of time for the case where the enhancement of the capture cross- 
section is incorporated. 



